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Abstract 

The numerical calculations of the standard Eliashberg-Nambu strong 
coupling theory and a formula of isotope effect derived in this paper 
provide direct evidences that an-harmonic vibrations of lattice en- 
hance isotope effect with anomalous coefficient a > 1/2. The results 
in this paper explain very well the wide distributed a values for the 
samples with different ratios of substitutions of ^^C by ^^C of fulleride 
superconductor RbsCgo- The calculations of isotope effects indicate 
that the intra-molecule radial modes have more important contribu- 
tions to superconductivity than the intra-molecule tangential modes 
with higher phonon frequencies 
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1 Introduction 

The discoveries of isotope effects in metallic superconductors confirm the 
electron-phonon mechanism as the origin of attractive interaction between 
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electrons. Anomalous isotope effects such as the negative isotope effects 
with a < and the enhanced isotope effects with a > 1/2 are found in many 
new types of superconductors. Understanding the origin of these anomalous 
isotope effects is very important to understand the microscopic mechanism 
of superconductivity. 

Coefficient of isotope effect 0<Q;<l/2is correct for most conventional 
superconductors. Negative isotope effect a < 0, where Tc increases with mass 
M, had been found experimentally in conventional superconductors such as 
PdH [Il[2] (a = —0.25), organic superconductors ^|4j and a-Uranium (a = 
—2.0) [5]. The negative isotope effect had also been found in Sr2Ru04 [6] 
when ^^O atoms were substituted by ^^O atoms. 

Large positive Oxygen isotope effects a > 1/2 beyond the BCS value 
1/2 had been found in HTSC material La2-xSr2;Cu04 (a=0.75) at doping 
level x=0.12 \7\ and in fuUerides RbsCeo (« > 1-0) with different ratios of 
substitutions of ^^C with [8l[9lll0|. For La2-xSr^Cu04, near x=0.12 the 
large coefficient of isotope effect is probably due to the strong an-harmonic 
vibrations companying with structural transitions. Similarly, the transition 
from partial substitution to completed substitution of ^^C by ^^C for RbaCeo 
leads to the decrease of a from 1.275 or 1.189 ^i9j to 0.30 [11] because the an- 
harmonic effect is relatively weak at completed substitution. The Eliashberg 
theory combined with different models of an-harmonic lattice vibrations had 
been used to qualitatively explain these anomalous isotope effects [121 [131 HI] . 
Beyond Migdal's theorem, the coefficient of isotope effect a > 1/2 had also 
been obtained by including non-adiabatic effect [15j. 

It's very convenient to find a formula for isotope effect which is able to 
explain anomalous isotope effects such as the enhanced isotope effects with 
a > 1/2 and large negative isotope effect. Most of previous formulas of iso- 
tope effects don't include the frequency shift due to the isotope substitution. 
A formula for isotope effect is derived in this paper using McMillan's Tc 
formula, which explicitly includes the frequency shift of phonon and is qual- 
itatively consistent with the more accurate numerical results of Eliashberg- 
Nambu theory. The formula is similar to a previous formula derived with 
similar idea [IB]. In this work, we not only compare the results of the formula 
with the strong-coupling theory, but also we study in details the meaning of 
the an-harmonic coefficient A* and its relation with the potential V{r) of 
the vibrations of atoms around their equilibrium positions in a crystal. We 
concentrate on the anomalous isotope effects of fuUeride superconductors [17] 
in the present paper. 
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An important relation is the key to understand isotope effect and written 

as 

M{uj^)X = r] = const (1) 

where the McMillan-Hopfield parameter r] characterizes the chemical envi- 
ronment of atoms in a material [HI [19] and is expressed as rj = N{0){J'^), 
where A^(0) is the density of state at Fermi energy, J the matrix element of 
electron-phonon interaction. Certainly, we assume t] is constant for isotope 
substitution or others small structural changes. We will see in section (jlj) 
that, under the constraint, the an-harmonic effect is equivalent to the fact 
that the parameter A is dependent on M the mass of atom. In a previous 
work, using electron-phonon mechanism and the constraint Eq.([T]), we have 
successfully explained the spatial anti-correlation of the energy gap and the 
phonon energy of Bi2212 superconductors observed in STM experiments pU] . 
In this work, using the constraint combining with McMillan's Tc formula, we 
derive a formula of isotope effect that can explain almost all anomalous iso- 
tope effects, especially for fulleride superconductors. 

2 Theory 

The energy-gap equation of Eliasgberg-Nambu theory in the Matsubara's 
imaginary energy form is standard [T5|. When temperature is very close to 
transition temperature Tc, the Eliashberg equation can be simplified as 

N 

J2{Kmn- p5mn)^n = (m > 0) (2) 
n=0 

K^ri = A(m-n) + A(m + n + 1) -2/i*(A^) 

m 

- 5^„(2m + l + A(0) + 2 5]A(/)), 

1=1 

where A„ = An/|a;„| is the energy-gap parameter and X{n) = 2 /q°° dLja'^F{uj)uj / {00"^+ 
(27rnT)^). The pair-breaking parameter p is introduced to form a eigenvalue 
problem with p=0 corresponding the physical energy-gap equation. The tran- 
sition temperature Tc is defined as the temperature when the maximum of 
eigenvalue of kernel matrix Kmn crosses to zero and changes its sign. We use 
about N=200 Matsubara energies to solve above equation. 
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The Eliashberg function is expressed as 



( 



c 







(a;3)2+{L^2)2' 



c 



UJ — ^lp\ < UJ3 

others, 



(3) 



where f2p is the energy (or frequency) of phonon mode, UJ2 the half-width 
of peak of phonon mode and u)^ = 2uj2. We can write the parameter of 
electron-phonon interaction A = A(0) = 2 duja^F{uj) / uj. The moments 
(a;") of the distribution function {2/ X)a'^F{ijj) /oj are defined as (a;") = 
2/A /q°° duja'^F{uj)uj'"'~^ . The parameter A characterizes the strength of electron- 
phonon interaction by A oc N [0) {J"^) / M {uj"^) in terms of Eq. ([1]) . The Coulomb 
pseudo-potential is defined as /xq = N{0)U and its renormalized value as 
fi* = /io/(l + /^o lii(-^c/i^o))) where U is the Coulomb parameter, Eq the 
characteristic energy for electrons such as the Fermi energy or band width, 
and Uq the characteristic phonon energy such as the energy cutoff of phonon 
energy or Debye energy. If c<j2 ^ ^p, Eq.([T]) can be simplified as 



which is obviously correct for Einstein model a'^F{uj) = {X/2)Qp6{u! — Qp). 

3 A formula of coefficient of isotope effect 

We will derive a formula of isotope effect by including an-harmonic effect. 
The shift of phonon energy (or frequency) is explicitly included in the for- 
mula. For single-frequency like mode, the a values in terms of the formula are 
qualitatively consistent with the numerical solutions of Eliashberg equation. 
Most importantly, the formula clearly shows that the an-harmonic vibration 
of lattice leads to the enhanced isotope effect with a > 1/2. It's conveniently 
to define a parameter A* hj Qp oc AI^^^^^*^/"^ to measure the an-harmonic 
effect. We can easily see that A* ^ represents the an-harmonic effect. The 
an-harmonic parameter A* can be expressed as 



Under harmonic approximation ^Ip (x M the parameter A is a constant 
{6X/6M = 0) because 



MfipA = 1] = const 



(4) 



A* = 1 + 2M/VLp{5np/5M). 



(5) 



A* = 1 + 2M/VLp{5VLp/5M) 



M/\{6\/6M) = 
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by using the relation Mf2pA = const. 

We start from the McMillan formula of transition temperature of su- 
perconductor 

Qd . 1.04(1 + A) ^ 
^^=L45^"Pf-A-/.*(l + 0.62A)]- 

where A is the strength of electron-phonon interaction. The coefficient a of 
isotope effect defined by Tc oc can be obtained from above McMillan 

formula by the direct mass-dependent from Gd oc M~^/^ and implicit mass- 
dependent from fi* by coq oc M~^/^. If an-harmonic effect is included then 
Qj:,,ujq oc M"*^^"^*)/^. In this work we consider additionally mass-dependent 
from A by the well known constraint MfipA = r], we obtain a formula of 
isotope effect which is expressed as 

1 _ 1.04(1 H-A)(l+0.62y 

2 2[A + 0.62A)]2 ^ ''^ ^ ^ ^ 

where A* is defined in Eq.(j5]) and the function T(A,/i*) is expressed as 

1.04A(l + 0.38/i*) 
[A + 0.62A)P 



or 



_ A(2.08-A) + A(2.7904 + 1.24AK 
^^^'''^ - 2[A-^*(1 + 0.62A)P 
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(0.04 + 0.42A)(l + 0.62A)(/i*; 
^ 2[A-/i*(l + 0.62A)]2 ' 

dependent on whether Qd^uq oc 

for Eq.® or Bd, 

for Eq.(l9]). We notice that, for the Eq.Q, the an-harmonic effect enters 
into the coefficient a only by the M-dependent A, however for the Eq.Q not 
only by M-dependent A but also by the 6d and loq. It's very important that 
T{X,fi*) > for reasonable values < A < 2 < 2.08, so the sign of the 
third term of Eq.([7]) is determined by the an-harmonic parameter A*. It's 
obviously that A* = under harmonic approximation f2p oc M~^/^ and the 
third term in the Eq.Q is also zero. The derivative 6Qp/6M used in the 
calculation of an-harmonic parameter A* is approximately obtained from the 
experimental energy (or frequency) shift of phonon. The first two terms give 
a < 1/2. The properties of the third term to isotope effect is determined 
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by the sign of the an-harmonic parameter A*. Generally, Qp decreases with 
increasing mass M, thus 2{M /VLp)5VLp / 5M < 0. If 2\{M /VLp)5VLp / 5M\ < 1, 
A* > 0. So we can get a > 1/2 if the third term has larger absolute value 
than the second term. 

4 The anomalous isotope effects with a > 1/2 

In many literatures of A„C6o, the average mass m of Ceo is used in the 
definition of coefficient of isotope effect by d = —{m/Tc)6Tc/6m. If the 
ratio of substitution is p, the average mass of Cqo molecule m = 60[Mi2(l — 
p) + Misp], so Sm = m- 6OM12 = 60p5M with SM = M13 - M12. The 
relation a = a/p connects the parameter a with the usual definition of 
a = —{M/Tc){6Tc/6M). The usual coefficient of isotope effect a = pa = 
1.275 is obtained with a=2.125 and p=0.60 [8], a = 1.189 with a=1.45 
and p=0.82 [9], a = 0.462 with a=1.4 and j9=0.33 [10]. For completed 
substitution with p=1.0, a = a=0.30 [TT] . 

The relation M{uj'^)X = 77 = const services as a constraint to determine 
the parameters A and Qp in numerical calculations. The Coulomb pseudo- 
potential /i* has to change in isotope substitution because /i* is dependent 
on the cutoff of phonon energy caJq oc M~^/^ or oc M"^^""^*)/^ in an-harmonic 
approximation. In harmonic approximation, 6fi* = — (/i*)^5M/2M, the an- 
harmonic effect is included in calculations only by the M-dependent A. The 
an-harmonic effect can be realized by shifting phonon energy 6Qp to make 
A* = 1 + 2M/VLp{8VLp/5M) ^ 0. 

The intra- molecule radial mode, which is about f2p=65 meV or 525 
cm~^ from infrared spectrum [inilB] and the intra-molecule tangential modes 
around 1400"^ or 174 meV |T7j have strong intensity. However we concentrate 
attentions on the intra-molecule radial mode with energy VLp=Qb meV and 
half-width u;2=8 meV. The Coulomb parameter /io = f/A^(0) and the corre- 
sponding renormalized Coulomb parameter yU* = /io/[H-/io ln(-Ec/i^o)] can be 
estimated from ab-initio density functional theory (DFT) based on pseudo- 
potential method using atomic orbital basis functions [21] • In the DFT cal- 
culation, the super-cell includes one Cqo molecule and three Rb atoms. The 
effect of orientation of different Ceo molecules is ignored. The parameter U is 
the charge energy defined as S'^Etotl^in? = E{n + 1) + E{n — 1) — 2E{n), E{n) 
the total energy of electric- neutral system of n valence electrons, E{n + 1) 
and E{n — 1) the total energies with one negative and one positive charge 
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Table 1: The Coulomb parameter /i* is calculated when we choose valence- 
electron configuration (a) Ap^5s^ for Rb atom and 2s^2p^for Carbon atom 
in DFT calculations. The case (b) is corresponding to the general value 
/x*=0.10. The coefficients a of isotope effects are calculated using the 
phonon-energy shifts from 65 meV to 62.3 meV after the isotope substitution 
(A*=0.003). 





U(eV) 


N(0)(l/eV) 


UN(0) 






A 


Tc(K) 


a 


(a) 


1.991 


2.197 


4.374 


6.9078 


0.127 


0.700 


29.544 


0.342 


(b) 










0.100 


0.670 


29.820 


0.358 



respectively. We choose the possible valence-electron configuration 4p^5s^ 
for Rb atoms and 2s^2p^for Carbon atoms. The core electrons are presented 
by Troulier-Martins pseudo-potentials. The electrons in semicore state 4p 
of Rb atoms having already treated as valence electrons have the single-^ 
basis set and all others valence electrons for all atoms have the split valence 
double-^ plus polarized basis sets. We use T point sampling the first Brillouin 
zone. The exchange-correlation potential is GGA Perdew-Burke-Ernzerhof 
type [22] and the spin-polarization effects are included in the self-consistent 
calculations. The results in table([T]) show that the Coulomb parameter hq is 
equal to 4.374 and the renormalized Coulomb parameters fi* is 0.127 if the 
maximum of phonon energy ujq is 150meV and Ec=15 eV when all valence 
electrons are included. The common value /i*=0.10 is close 0.127 obtained 
in this work. We have known the Coulomb parameter n* and the phonon 
energy Qp, the parameter A of electron-phonon interaction can be defined by 
the experimental transition temperature. From table([T]), we can see that the 
parameters fi* have small infiuence on isotope effects. Below we present the 
calculations of coefficient a of isotope effects with /i*=0.1. 

The parameter of electron-phonon interaction A=0.67 is defined by the ex- 
perimental transition temperature T^ about 29.5(K) using parameters /i*=0.1 
and f2p=65 meV. We calculate the t] = rii2 for ^^C. The phonon energies after 
^^C substitutions are unknown because the phonon energies are dependent 
on the ratios of ^^C substitutions. We choose the nine possible energies 
for phonons in table ([2]), which are all smaller than flp=Q5 meV because 
the energies (or frequency) decrease with increasing mass of atoms. We 
also assume that the parameters A of electron-phonon interaction alter af- 
ter isotope substitutions. The new Eliashberg function a^F{uj) after isotope 
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Table 2: The simulation parameters and the results with different an- 
harmonic parameters A*. The coefficients a are calculated based on the 
numerical solutions of Eliashberg equation and a' from the formula Eq.([7]) 



and Eq.( 


ID 


derived in this 


paper. 










i^C M (u) Qp (meV) 


A 


TciK) 


A* 


a 


a' 


(1) 


13 


64.3 


0.632 


25.5 


0.742 


1.929±0.10 


2.377 


(2) 


13 


63.4 


0.650 


27.2 


0.409 


1.207±0.09 


1.516 


(3) 


13 


62.3 


0.674 


29.0 


0.003 


0.358±0.09 


0.464 


(4) 


13 


62.2 


0.676 


29.2 


-0.034 


0.282±0.09 


0.368 


(5) 


13 


62.0 


0.680 


29.5 


-0.108 


0.132±0.08 


0.177 


(6) 


13 


61.8 


0.685 


29.8 


-0.182 


-0.02±0.08 


-0.014 


(7) 


13 


61.6 


0.689 


30.3 


-0.255 


-0.165±0.08 


-0.205 


(8) 


13 


61.2 


0.698 


30.9 


-0.403 


-0.457±0.08 


-0.588 


(9) 


13 


60.0 


0.727 


33.2 


-0.846 


-1.305±0.08 


-1.736 


12C 


12 


65.0 


0.670 


29.82 









substitution by ^^C is obtained from the old one before isotope substitution 
by simply shifting energy of phonon and scaling it to satisfy the constraint 
M{u^)X = ?7i3 = ?7i2. Based on the new Eliashberg function a^F{uj), we 
can calculate the Tc and A after the isotope substitution. The new Coulomb 
pseudo-potential after isotope substitution /i*=0. 099583 is obtained from the 
formula = — (yLi*)^5M/2M. The nine possible values of A corresponding 
to nine possible energies of phonon are collected in table ([2]) . The transition 
temperatures Tc are obtained by solving Eliashberg equation Eq.([2]). The 
coefficients a are easily calculated in terms of the transition temperatures Tc 
after and before isotope substitutions by a = — ln(^^Tc/^^Tc)/ln(Mi3/Mi2). 
The an-harmonic parameters A* are calculated in terms of the shifts of 
phonon energies. 

For the case (2) in table ([2]), we get a = 1.207 which is very close to 
the already known maximum value 1.275 in experiments with uncompleted 
substitution (a=2-2.25, p=0.60) [Sj. For sample with completed substitution 
with ^^C, the value of a decreases to ~0.30 [IT]. To explain the interesting 
results, we assume the increase of ratio of substitution makes the distribution 
of ^^C more homogenous and the an-harmonic effect become weak with small 
A*. If the phonon energy shifts to 62.3 meV with smaller an-harmonic pa- 
rameter yl*=0.003, the value of a=0.358 is close to the experimental values 



8 




Figure 1: The comparison of coefficients of isotope effect having obtained 
from numerical calculations of Eliashberg equation and from Eq. ([7]) for radial 
intra- molecule mode Qp=Q5 meV (a) and the tangential intra- molecule mode 
np=150 meV (b) 



from 0.30 to 0.37 for the samples with completed substitution. 

For negative an-harmonic parameter A* =-0.1 82, we find a=0.02 close to 
zero. After having made A* more negative further, at 74*=-0.403, we get the 
negative coefficient a=-0.457. However, the negative isotope effect wasn't 
found in fuUeride superconductors. From above calculations we can see that 
the sign of A* determines the sign of a at relative larger absolute values 
of A*. If the absolute values of an-harmonic parameters are not too large 
0.07 > A* > —0.21 or an-harmonic effect is weak, the coefficient a are 
from to 0.5 within the range of general values. The above calculations 
explain the wide distributed values of coefficients for RbsCgo in experiments. 
This is because the an-harmonic parameters are dependent on the ratios of 
substitutions. 

From table ([2]), we can see that the values obtained from numerical cal- 
culations are very close to the values from Eq.d?]) with Eq.®. The values 
of a and a' have the same sign on the same rows. It's very important that 
if a > 1/2, A* > must be satisfied. Thus, the enhanced coefficients a 
larger than 1/2 for RbaCeo in uncompleted substitution samples mean that 
the strong an-harmonic effects. If the an-harmonic effect becomes weak as 
approaching the completed substitution, the value of a will decrease and 
reach to smaller value about 0.358 at A* = 0.003. 
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^o^(2) 



(3) 



Figure 2: An-harmonic coefficients A* for the models of interaction with 
form V{r) = hr'^. The inserted figures (1),(2) and (3) illustrate the shapes of 
interacting potentials and the quantized energy levels for h=0.01 and d=0.5, 
2.0 and 4.0 respectively. The maximum of quantum number of angular mo- 
mentum is up to /max =8 



Fig. [11(a) shows the results of the numerical solutions of Eliashberg equa- 
tion and the formula Eq. ([7]) . We can see that Eq. ([7]) with Eq. ([8]) is more close 
to the numerical solution than using the Eq.Q although for Eq.Q more 
completed an-harmonic effects are included. We have preformed the same 
calculations for intra-molecule tangential mode with energy f2p=150 meV, 
/i*=0.285, A=0.67 and Tc=29.5(K). From Fig{T](b), there is larger slope of 
a-6Qp/Qp curve compared with the intra-molecule radial mode with energy 
Qp=65 meV. The coefficient of isotope effect of the intra-molecular radial 
mode q;=0.358 at A* ~0.0 is more close to experiments at completed substi- 
tutions a ~ 0.30 — 0.37. So it is more correlated with superconductivity of 
fuUeride than the intra-molecule tangential mode. 

5 An-harmonic parameter A* and models of 
lattice vibrations 

In this paper, we introduce an important parameter, the an-harmonic param- 
eter A*. To make parameter A* more realistic, we study the models of lattice 
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vibrations in a crystal. The atoms in crystal move around their equilibrium 
positions. The change of potential energy of an atom away from its equi- 
librium position is written as V{r) = hr'^ with radially local displacement 
r. The harmonic vibration is corresponding d = 2. We had numerically 
solved the Schrodinger equation [— ^^V^/2M + V{r)]ijj{r) = Etjj{r) under 
isotropic approximation, especially, the one- dimensional radial part converts 
into difference equation. We calculated the relation Eq{M) between energy 
of ground state and mass of atom for a series of discrete values of M. It's 
well known that for harmonic approximation Eq{M) oc 1/M°'^. The effec- 
tive an-harmonic parameters A* for interaction with the form V{r) = hr'^ 
are obtained by fitting the discrete Eq{M) functions with function Eq{M) oc 
l/M^^"^*)/^. From the Fig 121 if < 2 we can obtain A* > Q. So the enhanced 
isotope effects with a > 0.5 are hopefully found based on the formula Eq.Q 
obtained in this work. If d > 2 so A* < 0, we can get the normal isotope 
effects with a < 0.5 and the negative isotope effects with a < 0. We can 
see that the simple model V{r) = hr'^ is suitable for the isotope effects of 
fuUeride Rb3C6o for different ratios of substitutions because the an-harmonic 
parameters are within the range to obtain experimental values of a. 



6 Discussion and Summary 

The Figure 1(b) shows clearly that the formula Eq.([7]) isn't good approxima- 
tion to more accurate numerical solution for high-energy mode. The reason is 
probably that McMillan's formula isn't correct in some regions of parameter 
space. The McMillan's formula is good approximation when the parameter 
A of electron-phonon interaction is not too large (A < 1) and the phonon 
energy f2p is not too high. The inter-molecule phonon modes such as the 
vibration between Ceo molecules and the between alkali-metal atoms and Ceo 
molecules are ignored in this work because their energies generally smaller 
than 12 meV. To obtain Tc=30K, the parameters A of electron-phonon inter- 
action are at least 3.0 so the instability of lattice will destroy superconduc- 
tivity. However, these inter-molecule modes still have significant influence on 
the properties of fuUeride superconductors such as the differences of isotope 
effects for different substitution configurations i?63[(^^Ceo)x(^^Ceo)i-z] and 

In summary, the coefficients a of isotope effects obtained in this paper are 
very close to the values in experiments or within the range of experimental 
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values. The reductions of a with increasing the substitution ratios of ^^C 
are due to the reductions of an-harmonic effects of lattice vibrations when 
the ratios of substitution tend to 100%. The enhanced coefficients of isotope 
effects with a > 1/2 generally happen at the intermediate stage of transition 
from one phase to another. Finally, the formula Eq.([7]) and the numerical 
methods used in this paper are also suitable to study isotope effects of other 
superconductors if the electron-phonon interaction is the pairing mechanism 
for electrons. 
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